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ABSTRACT 

We present a numerical study of the evolution of molecular clouds, from their formation 
by converging flows in the warm ISM, to their destruction by the ionizing feedback of 
the massive stars they form. We improve with respect to our previous simulations by 
including a different stellar-particle formation algorithm, which allows them to have 
masses corresponding to single stars rather than to small clusters, and with a mass 
distribution following a near-Salpeter stellar IMF. We thus refer to the stellar particles 
as "stars" . We also employ a simplified radiative-transfer algorithm that allows the 
stars to feed back on the medium at a rate which depends on their mass and the local 
density. Our results are as follows: a) As in previous studies, the clouds formed begin to 
gravitationally contract long (several Myr) before star formation begins, but in turn, 
star formation begins several Myr before the largest scales terminate their collapse, 
b) Contrary to the results from our previous study, where all stellar particles injected 
energy at a rate corresponding to a star of ~ 10 M@, the dense gas is now completely 
evacuated from 10-pc regions around the stars within 10-20 Myr, suggesting that 
this feat is accomplished essentially by the most massive stars, c) To further test this 
suggestion, we performed two additional simulations where the IMF "saturates" at 8 
and 20 M@. In the former, the clouds are never dispersed, while in the latter they 
are, although over a larger timescale than in the case with a full IMF. d) At the scale 
of the whole numerical simulations, the dense gas mass is reduced by up to an order 
of magnitude, although star formation (SF) never shuts off completely, indicating 
that the feedback terminates SF locally, but triggers new SF events elesewhere, over 
distance scales of a few tens of pc. e) The SF efficiency (SFE) is maintained at the 
~ 10% level, depending on the geometry, even though the gas continues to always 
infall onto the local star forming regions. We speculate that the inclusion of supernova 
explosions and magnetic fields should bring the SFE down to the few-percent level. 

Key words: interstellar matter - stars: formation - turbulence 



1 INTRODUCTION 

Understanding how star formation (SF) proceeds in our 
Galaxy, and other galaxies in general, is a key quest in as- 
trophysics. In recent years, it has become clear that the 
evolution and distribution of SF in the Galaxy is inti- 
mately linked to the structure and evolution of the molecular 
clouds (MCs) where it tak e s place (see, e.g., the reviews by 
Mac Low fc Klessenl 120041: iBallesteros-Paredes et all 120071 : 



l Audit fc Hennebelie] 
Audit & Hennebelie 



2005; 
2010; 



iHennebelle fc Audit] 
Hennebelie et al. 



2007; 



2008 



Hcits ch et alj|20 05; Heitsch et al 2006; Hcits ch fc Hartmannl 
l200Sl : iHeitsch et al.l 120091 ; iBaneriee et al.l 120091 ). the evo- 
lution of MCs and of their SF activity, from their 
formation by condensation of the atomic gas in the in- 
terstellar medium to their star-forming stages and, in 
some studies, to their destruction by stellar feedback 



Mac Low &z Klcsscn 2004: tsallcstcros-rarcdcs ct al. 2007: try, a — j — : — 7 — ninr.irJ\ a ,, , ,■ 

„ , T , „ ^ ., in' | p' ' „ , J|„„,J — , f Vazquez-Scmadcm ct al. 2010 . Among the above studies, 

McKeefcOstrikiJia |Va Z que Z -SemadelJ[20ia, and ref- \^ ,,,,,_ „„„?!, KftZZZ q11^ ; ZL .u ll onnJ 



erences therein). 

In a ser ie s of p r evious papers (|Vazquez-Semadeni et al.l 

l200d . 12001 l20ld . |2011|1. we have invest i gated along- 
side other groups (|Hennebelle fc Peraurd Il999l . |2000| ; 



those including self-gravity (IVazquez-Semadeni et al 



20081 ; IHeitsch et alj|200gf ) 



I l2007l . 

l2ooah 



l20ld . l201ll ; IHeitsch fc Hartmann 
have shown that the coherent formation of large clouds 
(several tens of parsecs) leads to the onset of global gravi- 
tational contraction throughout the cloud, at a stage when 
the cloud is still mostly composed of atomic hydrogen, with 



* e-mail: p.colin@crya.unam.mx 
© 0000 RAS 



Colin et al. 



SF only starting several Myr later, when the cloud has 
become mostly molecular. 

This result, however, is in contradiction with the largely 
established notion that MCs cannot be collapsing freely, 
since otherwise their resulting SF rates (SFRs) would 
be up to two orders of m agnitude larger than observed 
IjZuckerman fc Palmerlll974r ). A related property is that the 
observed SF efficiency (SFE, the fraction of a cloud's mass 
that ends up in stars) for whole giant molecular clouds 
(GMCs) is estim ated to be of only a few percent (e.g. 
iMvers et al.ll 19861 ). Hence, it is generally believed that MCs 
must be in or near equilibrium, supported against their 
self-gravity by supersonic turbulence, magnetic fields, or 
some combination thereof. Specifically, a number of SF the- 
ories have appeared in recent years in which the underly- 
ing scenario is that MCs are supported globally by turbu- 
lent pressure, while gravitational collapses occur only lo- 
cally, caused by the supersonic turbulent compressions (e.g. , 
Padoan fc Nordlundl200ll201ll ; lilrumholz fc McKedl2005l ; 
Hcnncbelle & Chabricr 2011, see also the discussion by Fed- 
errath & Klessen 2013). 

Of course, an alternativ e explanation has been known 
for over four decades (e.g ., Field 19701; Whitworthl 19791; 
Elmegreen|[l983l : ICoxlll983l : iFranco. Shore fc Tenorio-Taglei 
1994 ) for the low observed SFE of GMCs, namely that stellar 
feedback, mainly from the ionizing radiation from massive 
stars, may disrupt the clouds before they have converted 
much of their mass into stars. In this scenario, there is no 
need to support the clouds against their self-gravity Also, 
this scenario becomes even more feasible in view of the re- 
cent result, from numerical simulations of cloud formation 
and evolution, that the collapsin g clouds undergo hierar- 
chica l gravitational fragmentation (IVazquez-Semadeni et al.l 
2009). That is, the turbulent density fluctuations, having 
larger mean densities than that of the whole parent CMC, 
have shorter free-fall times and smaller Jeans masses, and 
therefore the densest clumps begin to form stars earlier than 
the rest of the cloud, and before the global collapse of the 
cloud t erminates. 

In IVazquez-Semadeni et al.l (|2O10l . hereafter Paper I), 
we carried out a first attempt to numerically capture this 
phenomenology, by performing simulations of cloud for- 
mation and evolution in the presence of ionization heat- 
ing from massive stars, which is considered to be the 
main feedback mechanism affecting GMCs of masses up 
to ~ 10 5 M G (iMatzner 1 12002| ; iKrumholz fc Matznerl 120091 ; 
iDale et al.l |2012j). In Paper I, the instantaneous, time- 
dependent SFE was defined as 



SFE(i) = 



M„{t) 



M dcnBO {t) + M,(t) : 



(1) 



where Mdcnse is the mass in dense gas (n > 100 cm -3 ) and 
M» is the total mass in stars. It was found in that paper 
that the prescription for ionization-heating feedback used 
there was able to maintain the SFE at the few-percent level 
throughout the evolution of the cloud, while control simu- 
lations not including it reached SFEs roughly an order of 
magnitude larger. Also, an analytica l model representing 
this sc enario was recently presented bv lZamora-Aviles et al.l 
(|2012f ). and shown to correctly describe several evolutionary 
properties of GMCs and their SF activity. 

However, one shortcoming of the feedback prescription 



used in Paper I was that it assumed that the stars responsi- 
ble for the feedback all injected energy into the medium at a 
rate roughly corresponding to that of a ~ 10 Mq star. This 
implied that the stellar feedback was possibly overestimated 
for clouds forming low-mass stars, and underestimated for 
clouds forming high-mass stars. In particular, Paper I found 
that the GMC-li ke clouds c o uld n ot be destroyed by the 
feedback. Instead. IDale et al.l i|2012l ). for example, have been 
able to disrupt clouds up to 10 Mq by means of ionization 
feedback. In Paper I, the SFE was kept low because the con- 
version of dense gas into stars was locally inhibited by the 
feedback, but not because the clouds at large were destroyed. 
Only the local clumps were destroyed. 

In the present paper, we improve on the numerical pre- 
scription used in Paper I in two ways. First, we use a prob- 
abilistic SF prescription instead of a fully deterministic one. 
As it turns out, this probabilistic prescription allows us to 
produce a mass spectrum for stellar particles, which can be 
tuned to resemble the Salpeter initial mass function (IMF). 
Second, once armed with a realistic stellar mass spectrum, 
we incorporate a mass-dependent ionization heating pre- 
scription for the feedback from the stellar particles produced 
in the simulations, applying a simplified description of ra- 
diative transfer, neglected in Paper I. We expect that, with 
this prescription, we can obtain a more realistic description 
of the effect of ionization feedback on MCs of various masses. 

The plan of the paper is as follows. In Sec. [2] we describe 
the numerical model, focusing in particular in the heat- 
ing and cooling functions employed (Sec. [2~2]l . the stellar- 
particle formation prescription (Sec. I2.3[l , the feedback pre- 
scription (Sec. 12. 4)1 . and the refinement criterion used (Sec. 
I2.1|l . We next describe the details and parameters of the 
simulations in Sec. [3] and then present our results in Sec. 
[4] In Sec. [5] we compare our results with those of Paper I 
and discuss some of their implications. Finally, in Sec. |6]wc 
present a summary and some conclusions. 



2 THE NUMERICAL MODEL 

The numerical simulations used in this work were per- 
formed using the adaptive mesh refinemen t (AMR) + INT- 
body Adap tive Refinement Tree code ART (jKravtsov et al.l 
1997, 2003). In the following sections we describe the adapta- 
tions we have performed to it for application to our problem 
of interest. 



2.1 Refinement 

The numerical box is initially covered by a grid of 128 3 (ze- 
roth level) cells. The mesh is subsequently refined as the 
matter distribution evolves. The maximum allowed refine- 
ment level was set to five, so that high-density regions have 
an effective resolution of 4096 3 cells, with a minimum cell 
size of 0.0625 pc. As in Paper I, cells are refined when the 
gas mass within the cell is greater than O.32M0. That is, 
the cell size is refined by a factor of 2 when the density in- 
creases by a factor of 8, so that, while refinement is active, 
the grid cell size Ax scales with density n as Aa; oc n -1 ' 3 . 
Once the maximum refinement level is reached, no further 
refinement is performed, and the cell's mass can reach much 
larger values. 
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Note that this constant-cell-mass refinement crite- 
rion does not confo rm to the so-called Jeans criterion 
ijTruelove et al.lll997l) of resolving the Jeans length with at 
least 4 grid cells. iTruelove et al.l i|1997l ) cautioned that fail- 
ure to do this might result in spurious, numerical fragmen- 
tation. However, we do not consider this a cause for concern 
since, as will be described in Sec. 12.31 our star formation 
prescription allows us to c hoose the s tellar- particle mass dis- 
tribution, and tune it to a lSalpeterl i|1955T ) value. 



2.2 Heating and cooling 

The main additional physical processes implemented in our 
simulations, and relevant to the physical problem studied 
here are a) the cooling and heating of the gas; b) its con- 
version into stars; c) the stellar feedback via ionization-like 
heating, and d) the self-gravity from gas and stars. 

We use heating (r) and cooling (A) functions of the 
form 



r 

A(T) 



2.0 x 10~ 26 erg s" 1 



10 exp 



+ 1.4 x 10" 



-1.184 x 10- 
T + 1000 



VTexp(^) 



(2) 



(3) 



These functions are fits to t he various heating and coo l- 



ing processes considered by iKovama fc Inutsuka 



uka 
sukal 



2000), 



2002) 



as given by equation (4) of Koyama & I nutsuka 
As noted in [y azqucz-Se madeni et al.1 (]2007 ). eq. (4) in 
IKovama fc Inutsukal (J2002J) contains two typographical er- 
rors. The form used here incorporates the necesary correc- 
tions, kindly provided by H. Koyama (2007, private commu- 
nication). With these heating and cooling functions, the gas 
is thermally unstable in the density range 1 < n < 10 cm -3 . 



2.3 Star formation prescription: a probabilistic 
approach 

In our simulations, SF is modeled as taking place in the 
densest regions, defined by n > nsF, where n is the gas den- 
sity, and «sf is a density threshold. If a grid cell meets this 
density criterion, then a stellar particle (SP) of mass m* may 
be placed in the cell, with probability P, every timestep of 
the coarsest grid. If the SP is created, it acquires half of the 
mass of its parent cell, and this mass is removed from the 
cell. Thereafter, the particle is treated as non-collisional and 
follows N-body dynamics. No other criteria are imposed. We 
set nsF = 9.2 x 10 4 cm -3 , which corresponds to a cell mass 
of 0.78 Mq at the highest refinement level. This fixes the 
minimum value for SP masses at 0.39 Mq. Note that, as in 
Pap er I, our SPs differ from the commonly- used sink parti- 
cles JBate et al.lll995l ; lFederrath et al.ll2010f l, mainly in that 
our SPs are not allowed to accrete after they form. Thus, 
we refrain from calling them as "sinks" , and use the nomen- 
clature "stellar particles" instead. However, as we describe 
below, our probabilistic approach to SP formation allows us 
to obtain a realistic mass distribution (IMF) for them. 

A few items are worth noting about our SF prescrip- 
tion. First, note that, once the maximum refinement level is 
reached, no further refinement is applied to a cell (cf. Sec. 
12. 1|) even if its density keeps increasing. Moreover, since the 
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Figure 1. Probability of having (solid line) and not having 
formed (dashed line) a stellar particle (SP) in a cell that meets the 
density criterion n > ngp after a certain number of coarse-grid 
time steps for P = 0.001. 



creation of an SP is a probabilistic event, the density of a 
cell where a gravitational collapse is going on continues to 
increase until an SP forms in the cell. Some authors have 
advocated the prescription that, once the maximum refine- 
ment level has been reached, a sink particle is created at the 
cell densit y that would correspo nd to the next refinement 
level (e.g., iFederrath et alJl201Cf ) in order to always fulfill 
the Jeans c riterion and thus com pletely avoid spurious frag- 
mentation l|Truelove et alJll997l ) until the sink particles are 
formed. However, we forgo of this recommendation since, 
as will be seen in what follows, our prescription allows us 
to impose the desired IMF of the SPs, and thus artificial 
fragmentation is not a concern. 

The prescription we use implies that the longer it takes 
to form an SP in a collapsing cell, the more massive the 
SP will be, because the cell's density will be higher. The 
probability of not having formed an SP after n s t ep s time 
steps is Vuo — (1 — P)" st °p B , while the probability of having 
formed it is V ycs — 1 - (1 - P) n = to P s . These probabilities 
are shown in Fig. [T]for P = 0.001. Serendipitously, we have 
found that the resulting mass distribution of the SPs is a 
power law, with an exponent that depends on the value of P. 
Thus, P is a control parameter that allows us to generate a 
stellar mass spectrum with the desired slope. In Figure[2]we 
show the evolution of the stellar mass spectrum in simulation 
LAF1 5; (cf. Sec.O with P = 0.003. The slope of the function 
changes from —1.21 at t — 25.6 Myr to —1.34 at the end of 
the evolution, thus hovering close to the Salpeter value of 
— 1.35. The most massive SP formed in this simulation has 
m* = 61 Mq, while the least massive ones have masses 
m* ~ 0.5 Mq. Note that we have no turnover of the IMF at 
small masses, but this is inconsequential for our purposes, 
since we are only interested in the feedback exerted by the 
stars on their parent cloud, and the low-mass stars exert no 
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Figure 2. Evolution of the spectrum of stellar masses for simu- 
lation LAFls;. The spectrum at any epoch can be well fitted by 
a power law. In the lower right panel, with a line, we also show 
this fit and the value of its slope. 



sig nificant feedback at the G MC scale (see, e.g., the review 
bv lVazquez-Semadenill201ll . and references therein) 

Finally, we note that, because now the SPs form in cells 
whose density is typically much larger than the threshold 
value nsF, in the present paper we choose tisf = 9.2 x 
10 4 cm -3 , to allow for a sufficiently large number of SPs 
to form. This is significantly smaller than the value used in 
Paper I, where SPs formed always at a density very similar 
to usf- Moreover, we stress that, contrary to the situation in 
our previous papers, our SPs now have masses correspond- 
ing to individual stars rather than to small clusters, and so, 
in what follows, we shall indistinctly refer to them simply 
as "stars". 



2.4 Feedback prescription 

Another important difference of our new feedback prescrip- 
tion, compared to that in Paper f, is in the way we im- 
plement the ionization feedback by massive stars. In Paper 
I, SPs injected thermal energy only to the cell where they 
were located (hereafter, the "stellar cell"), at a rate high 
enough to produce a realistic HII regiorp. Instead, here we 
now model the birth and evolution of HII regions by assign- 
ing a temperature of 10 4 K to all cells whose distance d to 
the SP satisfies the condition 



d<R s 



/_3 S»_ 

V 47rcm LO 



1/3 



(4) 



where _R S is the IStromgrenl |l939j) radius, S* is the flux of 
ionizing photons produced by the star, a = 3.0 x 10 -13 cm 3 



1 Because the thermal energy was dumped only in the cell were 
the SP was formed, neighboring cells were heated by numerical 
conduction. 



s _1 is the recombination coefficient, and tilos is a character- 
istic particle number density along the line of sight between 
the stellar cell and the test cell, which we discuss below. If it 
turns out that R s is smaller than the size of the stellar cell, 
we simply set the temperature of this cell equal to 10 4 K, and 
no further calculation is done. On the other hand, if Rs is 
larger than the stellar cell's size, then it is necessary to deter- 
mine whether d < R s or not. In principle, this poses a radia- 
tive transfer problem since, as is well known, eq. ((4]) is valid 
only for the case when the medium between the stellar and 
the test cells has a uniform density. However, if the medium 
is not uniform, then the photoionization-recombination bal- 
ance must be computed along the line joining t he SP and the 
grid cell in question (see, e.g.. |Pale et al.1120071 ). a procedure 
that can be quite computationally expensive. 

As a zeroth-order approximation to solve this problem, 
we opt for choosing a value for riLOS that can be deemed 
representative of the typical density along the path from 
the SP to the grid cell. Specifically, we take the geometric 
mean of the densities at these two locations in the simula- 
tion, tilos = \/nsp ntcst- This approximation for i? s is, of 
course, crude, and will miss, for example, shadowing effects 
due to intervening dense clumps between the stellar and the 
test cells but, for the purpose of modeling the large-scale 
dynamics of the MC containing these cells, we consider it is 
sufficient. We jokingly refer to this scheme as a "poor man's 
radiative transfer" (PMRT) scheme. 

For the ionizing flux S*, which depends on the 
SP's mass, we use t he tabulated data provided by 
iDiaz-Miller et al.l (|l99ST ). Note that only SPs with a mass 
greater than 1.9 Mq inject any significant ionizing feedback 
into the ISM, implying that only the massive SPs influence 
the dynamics of the molecular clouds. Finally, note that we 
turn off the cooling for the cells whose temperature is set to 
10 4 K. Otherwise, very dense cells would radiate away their 
thermal energy very quickly. Their temperature is held at 
10 4 K for a time t B , which we assume depends on the star's 



mass m» as 



2 Myr 

222 Myr (^) 



if m, < 8M ( 



'©; 



if m» > 8M« 



©■ 



(5) 



For stars more ma ssive than 8Mq, thi s time is a fit to 
the stellar lifetimes bv lBressan etal.1 (1 19931 ). while for stars 
with masses lower than that, it represents the fact that the 
duration of the stellar-wind phase is ~ 2 Myr, roughly inde- 
pendently of mass. This also means that we are representing 
the effect of the winds and outflows of low mass stars by an 
ionization prescription. While this is clearly only an approx- 
imation, we do not expect it to have much impact on our 
calculations, since the main source of feedback energy at the 
level of GMCs i s the ionization feedback from massive stars 
JMatznerll2002i 'l. 

Our prescription for the ionization feedback for massive 
stars was tested by running simulations of a box of 32 pc on a 
side, without self-gravity, filled with gas at uniform temper- 
ature and density, of 42 K and 100 cm -3 , respectively. These 
simulations used a resolution of 256 3 cells, with the adaptive 
refinement switched off. A massive SP (m* = 27 Mq) was 
placed in the center of the box, and the system was allowed 
to evolve freely. Figure [3] shows the expansion of the result- 
ing HII region, which is comprised of those cells with tem- 
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Figure 3. Evolution of the HII region produced by a 27-Mq star 
in a box of 32 pc per side. The black solid line represents the 
analytic solution given by eq. JB), while the blue dashed line is 
the radius of the HII region in the simulation. Temperatures of 
target cells located inside the Stromgren sphere, centered on the 
SP, are set to 10 4 K. The cooling is switched off in these cells 
during the lifetime of the star. 



peratures greater than a few thousand degrees. Als o shown 
in th is figure is the well known analytical solution i|Spitzer| 



7c s t\*/ 7 



«.(*)= ft (l+Ig 



(6) 



where Ri is the initial Stromgren radius and c s is the sound 
speed. The numerical solution is seen to agree with the an- 
alytic one to within ~ 30%, an accuracy we consider suffi- 
cient, given our interest only in the large-scale evolution of 
the clouds. 



3 THE SIMULATIONS 

Our simulations use the same initial setup as the runs in 
Paper I, which represents the evolution of a region of 256 pc 
per side, initially filled with warm gas at a uniform density 
of no = 1 cm -3 and a temperature To = 5000 K, imply- 
ing an adiabatic sound speed c s = 7.4 km s" 1 (assuming 
a mean particle mass /i = 1.27). The full numerical box 
thus contains 5.25 x 10 5 Mq. In this medium, we make two 
streams collide with a speed Vint — 5.9 km s _1 each (cor- 
responding to a Mach number of 0.8 with respect to the 
unperturbed medium) along th e x-direction (see Figure 1 of 
I Vazquez-Semadeni et alj 120071 ). The streams have a radius 
of 64 pc and a length of 112 pc each, so that the total mass 
in the two inflows is 9.0 x 1O 4 M0. Note that the streams are 
completely contained within the box, so that the compres- 
sion they produce is a single event. There is no continuous 
flow through the boundaries, as we use periodic boundary 
conditions. 

On top of the inflow velocity we superpose a field of ini- 



tial low-amplitude turbulent velocity fluctuations, in order 
to trigger the instabilities in the compressed layer that will 
cause it to fragment and become t urbulent IjHeitsch et al.l 
120051 ; I Vazquez-Semadeni et alJbOQq ). As in Paper I, we cre- 
ate this initial velocity fluctuati on field with a new ver- 
sion o f th e spectral code use d in IVazquez-Semadeni et al.l 
l|1995r) and lPassot et al.l |l995j), modified to run in parallel in 
shared-memory architectures. The simulations are evolved 
for about 40 Myr. 

The collision nonlinearly triggers a transition 
to the cold phase, forming a turbulent, cold, dense 



cloud (IHennebelle Sz Perault 
20051; 



Heitsch et al.l 



2005 



Vazquez-Semadeni et al.l 120061 ) 



19991; lAudit fc HennebeJiel 
iHeitsch et~al l200d : 

consisting of 



a com 
plex network of sheets, filaments, an d clumps of cold gas 



embe d ded in a warm diffuse substrat e 
20051; IHennebelle fc Inutsukal 120061 ; 



Audit fc Hennebelk 



Hcnncbcllc & Audit 



20071 ). The com plex as a whole quickly engage s in gravita- 



tional collapse (|Vazouez-Semadeni et al.l 120071 ). Moreover, 
the local density fluctuations become unstable and collapse 
in a shorter time than the global time because they are 
embedded in a contracting med ium and thus have shorter 
collapse times IjToala et al.ll2013T ). Eventually, they proceed 
to forming stars, which then heat their environment, 
forming expanding "HII regions" that tend to disperse the 
clouds. 

Although our results are based essentially on two sim- 
ulations, a few more runs were performed with a twofold 
purpose: to compare the old prescription for feedback to the 
new one, and to assert the importance of the most mas- 
sive SPs in the disruption of the clouds. All simulations but 
one were run using the "large-amplitude" (LA) initial ve- 
locity fluctuations (u rmB ~ 1.7 km s _1 ) as opposed to the 
"small-amplitude" (SA) case (i> rms ~ 0.1 km s _1 ) (see Pa- 
per I). Clearly, more fragmentation and more complex cloud 
structures are expected in the LA runs, thus causing them 
to produce somewhat smaller clouds that resemble low- or 
intermediate-mass star-forming clouds. On the other hand, 
the SA run allows us to consider a case of very high coher- 
ence and unifor mity, which tends to form a h igh-mass star- 
forming region (| Vazquez-Semadeni et al.ll2009l ). and further- 
more resembles the conditions wc used in previous papers 
(Vazquez-Semadeni et al. 2007, 2011). 

The nomenclature for the runs introduced in Table 1 
continues to use the acronyms LAF or SAF used in Paper 
I, where F stands for feedback. A number 1 or after the 
letter F means feedback is on or off, respectively. The test 
simulations are denoted with "old", 8 or 20 after the word 
LAF. 



4 RESULTS 

4.1 LAF1 and SAF1 runs 

The new schemes for the probabilistic star formation recipe 
and for the ionization feedback by massive stars were used 
to run the SAF1 and LAF1 simulations. Figures [4] and [5] 
show, for the full simulated box, the evolution of the mass 
in dense gas (n > 100 cm" 3 ), the mass in stars, the mass 
in both components, and the star formation rate (SFR) and 
efficiency (SFE), in the SAF1 and LAF1 runs, respectively. 
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Table 1. Parameters of the Simulations 



Run »; rm s 

name [km s -1 ] 



SAF1 


0.1 


LAF1 


1.7 


LAFO 


1.7 


LAFold 


1.7 


LAF8 


1.7 


LAF20 


1.7 



Feedback 



New prescription, full IMF 

New prescription, full IMF 

Off 

Old prescription from Paper I 

New prescription, max stellar mass = 8 Mq 

New prescription, max stellar mass = 20 Mq 



In contrast to Paper I, we now define the instantaneous SFE 

as 



SFE(i) 



M«{t) 



M max (i) + M,(i)' 



(7) 



where M max is the maximum mass in dense gas reached 
between the start of the simulation and time t. This is be- 
cause in the present simulations the clouds are eventually 
dispersed, and thus the SFE evaluated with the instanta- 
neous dense gas mass approaches unity at late stages of evo- 
lution, but not because all the dense gas has been converted 
to stars, but rather because the remaining gas is evaporated 
by the stars. Our definition gives us instead an approximate 
measure of the net SFE, that is, the fraction of the dense 
gas mass ever present in a given volume and over a certain 
time interval that is converted into stars. 

In both runs, the evolution of the mass in dense gas is 
similar (see the top left panels of Figs. U and [SJ: first, dense 
gas starts to accumulate as the evolution proceeds until it 
reaches a maximum, then the effect of the radiation feedback 
is such that it overcomes the buildup of dense gas by grav- 
itational accretion. Also, in both runs, SPs continue to form 
after the maximum in the mass in dense gas is reached, al- 
beit at a significantly declining rate. Later in the evolution, 
in the LAF1 run, the mass in dense gas begins to increase 
again. This time the density does not reach the threshold for 
star formation and thus no new SPs are formed (see bottom 
left panels of Figs. [4] and [5} . 

In run SAF1, the largest star-forming region forms in 
the center of the box, due to the coherent collapse of the 
entire sheet-like cloud formed by the collision. As in Paper 
I, we refer to this region as "the Central Cloud". We enclose 
the cloud in a cylinder of radius 10 pc and length 20 pc 
located at the center of the box. Figure[S]shows the evolution 
of the mass in dense gas, the SFE, the mass in stars, and 
the SFR for the Central Cloud. This figure is similar to Fig. 
2] (or Fig. 0) except that the bottom right panel now shows 
the evolution of the SFR instead of the evolution of the total 
mass. Unlike what happens with the whole box, where some 
dense gas still remains by the end of the evolution, here we 
witness the complete dispersal of the cloud from this region. 
In addition, in the lower left panel we see that the mass in 
SPs also decreases by the end of the evolution. Because our 
SPs have no winds and do not explode as supernovae, this 
can can only mean that the stellar cluster, formed from the 
dense gas mass of the cloud, is being dispersed as well (that 
is, its constituent stars are leaving the cylinder that initially 
contained the cloud). 

In the case of run LAF1, the clouds identified as Cloud 
1 and Cloud 2 in Paper I are also used here to study the 



100 

10 20 30 40 
l(Myr) 




\_.-r 1000; 



20 25 30 35 40 
t(Myr) 




10 20 
t(Myr) 



Figure 4. Evolution of the mass in dense gas (top left panel), 
the star formation rate (SFR; dashed line) and efficiency (SFE; 
solid line) (top right panel), the mass in SPs (bottom left panel), 
and the total mass (stars plus dense gas, bottom right panel) in 
the whole box in run SAF1. As a result of the destruction of the 
dense gas the mass in SPs as well as the SFE reaches a maximum 
before the end of the evolution. 



cloud evolution and star formation activity. Again, cylin- 
ders of the same size as that used with the Central Cloud 
are placed around Clouds 1 and 2. As in Fig. [5] Figs. [7] and 
[5]also show the evolution of the mass in dense gas, the SFE, 
the mass in stars, and the SFR for Cloud 1 and Cloud 2, 
respectively. As with the Central Cloud in run SAF1, Cloud 
1 is also destroyed, and no stars are left inside the cylinder 
where the cloud initially was; that is, the stellar cluster as- 
sociated with the cloud is dispersed, in ~ 12 Myr. Cloud 2 is 
also completely destroyed, but unlike Cloud 1, here we still 
can find SPs inside its corresponding cylinder, as it was the 
case for the Central Cloud. Figure [5] shows the evolution in 
the neighborhood of Cloud 2 over 15 Myr, illustrating these 
results. 



4.2 Feedback scheme comparison 

In the feedback scheme of Paper I, HII regions were created 
by the injection of thermal energy from SPs. The energy 
was deposited entirely in the cell where the stars were lo- 
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Figure 9. Cross-section images of the density field in the neighborhood of Cloud 2 in run LAF1, at times (in Myr) 20.59 (top left), 
25.62 (top right), 30.1 (bottom left), and 34.86 (bottom right), showing the dispersal of the cloud. The black dots show the stellar 
particles (SPs). The horizontal ruler shows a scale of 26.3 pc. Note that the density field is shown on an inclined cross section through 
the simulation, but the SPs are shown in 3D space, so all particles in front to the density plane can be seen. Note the complete dispersal 
of the cloud within 15 Myr. 



cated, and thus neighboring cells were heated exclusively by 
conduction, rather than by radiative heating. The value of 
the rate at which the energy was dumped was chosen so 
as to produce reasonably realistic HII regions. Additionally, 
the cooling in the heated cell was turned off, since otherwise 
most of energy would be radiated away in these initially very 
dense cells. Thus, it is not feasible to directly compare the 
results of our new simulations, with the PMRT scheme used 
in the present paper. However, it is important to compare 
the old prescription with the new one used in the present pa- 
per in a controlled manner, to assess the differences induced 
by the prescription, in addition to the differences induced 
by the presence of a stellar IMF. Therefore, we have run an- 
other LAF-type simulation, labeled LAFold, which uses the 
old feedback prescription from Paper I, and in which all SPs 
with M* > 10 Mq inject thermal energy at a rate equal to 
that used in Paper I. 

Figure [TTJ1 shows the evolution of the mass in dense gas, 
the SFE, the stellar mass, and the total (dense gas plus stars) 
mass for runs LAF1 (black solid lines), LAFO (red dotted 
lines) and LAFold (blue, short-dashed lines), together with 
two other runs to be discussed in the next section (see Table 
[TJ. We see that the LAFold run does not disrupt the clouds, 
in line with the results of Paper I. 



4.3 The role of the most massive stars in the 
destruction of the clouds 

To assert the importance of the feedback of the massive stars 
in the destruction of the clouds, two extra LAF models were 
run, labeled LAF8 and LAF20. LAF8 (magenta long-dashed 
lines in Fig. I10p is a run with the new feedback prescription, 
but with all SPs with M* > 10 Mq ionizing their surround- 
ings as if they were a star of 8 Mq. LAF20 (cyan dot-dashed 
lines), moreover, is a run similar to LAF1 but in this case 
the feedback "saturates" at 20 Mq; that is, all SPs with 
M* < 20 Mq have the feedback they should according to 
their mass, but those SPs with M* > 20 Mq exert a feed- 
back as if they were a star of 20 Mq . 

From Fig. [TO] we see that, like run LAFold, run LAF8 
does not destroy the clouds; the mass in dense gas in the 
whole simulated box (and in the individual clouds, not 
shown) continues to increase. Run LAF20 is an intermediate 
case between those runs in which clouds are not destroyed 
and run LAF1: in run LAF20, the mass in dense gas reaches 
a peak before the end of the evolution. Because the feedback 
in run LAF20 is not as strong as it is in run LAF1, this max- 
imun is reached few Myr later. This experiment demostrates 
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Figure 5. Evolution of the mass in dense gas (top left panel), the 
SFR (dashed line) and SFE (solid line; top right panel), the mass 
in SPs (bottom left panel), and the total mass (stars plus dense 
gas, bottom right panel) in the whole box for run LAF1. Unlike 
the SA run, here the star formation does not stop abruptly, and 
instead a slow decline in the star formation rate is observed. 
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Figure 6. Evolution of the mass in dense gas (top left panel), 
the SFE (top right panel), the mass in SPs (bottom left panel), 
and the SFR, (bottom right panel) in the Central Cloud in run 
SAFl. Ioinization feedback is so efficient that the cloud only lives 
abut 10 Myr. Interestingly, the stellar cluster is almost dispersed 
in the next 10 Myr of evolution. 



that stars with M* > 20 Mq are crucial for the destruction 
of clouds of masses up to a few times 10 4 Mq. 



5 DISCUSSION 

5.1 Comparison with previous work 

The effect of feedback has been studied by numerous work- 
ers, bot h analytically and nume rically (see, e.g., the re- 
view by IVazauez-Semadenill201ll . and references therein). 
In particular, the pioneering numerical simulations of 
iBania fc Lvonl (1980) included various cases of heating and 
cooling functions for the medium, and used radiative trans- 
fer (on a 40 x 40 two-dimensional grid) to include the ef- 
fects of photoionization from OB stars in a 180-pc square 
region, making them a direct precursor of this work. Even 
at their very limited resolution, they foresaw several out- 
comes of this setup, such as the formation and maintenance 
of a cloud population, that the clouds would be gravitation- 
ally unstable had self-gravity been included, and that the 
SFR would be self-consistent if 0.1-0.5 of the mass in the 
clouds were to go into the formation of new massive stars, 
thus making a prediction for the SFE. However, their re- 
quired SFE for self-consistency was too high, presumably 
because self-gravity was not included in their simulations. 
This limitation also meant that the feedback stars had to 
be placed randomly in the simulation. Subsequent numer- 
ical works focused mostly on the effect of supernova feed- 
back on the structuring of the ISM on kiloparsec scales (e 



Rosen et al 



19931; iRosen fc Bregmanlll995l;lde Avilledfe OQO; 
Mac Low et al.l|2005l; iJoung fc Mac Lowil200tj ; IWood et all 
20101 ; iHill et alj|2012h but, in general, self-gravity has not 
been included in these works, and the supernova rate has 




25 30 

t(Myr) 



Figure 7. Evolution of the mass in dense gas (top left panel), 
the SFE (top right panel), the mass in SPs (bottom left panel), 
and the SFR (bottom right panel) in Cloud 1 of run LAF1. 



been an input parameter for the simulations, rather than a 
self-consistent output. 

Another line of study has been the simulation of 
feedback by stellar outflows at the clump (parsec) scale, 
aiming at either maintenance of the turbulence within 
the clumps (e.g.. iLi fc Nakamurall2006l ; I Carroll et all 120091 ; 
ICunningham et al.ll2009 | ). or at the self-regula tion of star 
formation (e.g.. | Li fc Nakamurall20 06; Naka mura fc Lill2007l ; 
IWang et al.l 120101 ). The latter are closest in aim to our 
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Figure 8. Evolution of the mass in dense gas (top left panel), 
the SFE (top right panel), the mass in SPs (bottom left panel), 
and the SFR (bottom right panel) in Cloud 2 of run LAF1. 
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Figure 10. Evolution of the mass in dense gas (top left panel), 
the SFE (top right panel), the mass in SPs (bottom left panel), 
and the total mass (stars plus dense gas, bottom right panel) in 
the whole numerical box of each of the five LAF runs in the whole 
box: run LAF1 (black solid lines); run LAF0, with no feedback 
(red dotted lines); run LAFold (blue short-dashed lines), with the 
feedback prescription from Paper I; run LAF8 (magenta long- 
dashed lines), with the new feedback prescription but with all 
SPs with Af« > 10 Mq ionizing their surroundings as if they had 
a mass of 8 Mq; and finally, run LAF20 (cyan dot-dashed lines), 
in which all SPs with Af* < 20 Mq feed back according to their 
masses, but SPs with M„ > 20 Mq feed back as if they were a 
star of 20 Mq . 



present study, although not in scale, as they only consider 
parsec-sized regions, and only outflow feedback, which cor- 
responds to the effect of low- and intermediate-mass stars, 
which does not seem to be the dominant driver at the GMC 
scale (|Matzner ||20o3 >. 

Our re sults in thi s pape r are most directly comparable 
to those of iDale et al.l IJ2012T ). who performed a parameter- 
space study of the disruptive effect of photoionizing radia- 
tion of molecular clouds of various masses. Thus, the basic 
physical processes at play in their simulations are very sim- 
ilar to those included in ours. Their main result is that, 
while clouds of masses up to ~ 10 5 Mq can be readily de- 
stroyed by the ionization feedback from their newly-formed 
stars, clouds with M ~ 10 6 Mq cannot be destroyed, as 
their escape velocities are larger than the sound speed in 
the photoionized gas. This regime is not sampled by our 
simulations, in which the total dense gas mass is never more 
massive than a few times 10 4 Mq. 

The main differences between our setup and theirs are 
that they use a polytropic equation of state covering only 
a temperature range corresponding to molecular and cold 
atomic gas, and that they start with a suite of initial clouds 
in various configurations, while we let the clouds form self- 
consistently out of the warm ISM. Also, they include the 
photoionized gas resulting from the stellar feedback, but 
their simulations lack the warm (neutral and ionized) sub- 
strate in which the clouds dwell, and out of which they form 
in our simulations. That is, in their simulations there is no 
possibility of the warm environment penetrating in to the 
clouds , as proposed theoretically by Hcnncbelle & Inutsuka 
(2006), and suggested observationally bv lKrco et al.l Q2008). 
Thus, our self-consistently formed clo uds may b e mor e 
porous, and thus less bound, than those of lDale et al.l (J2012J). 
Also, our probabilistic SF prescription allows our SPs to 
be individual stars always and with a realistic IMF. This 
means that our simulations can be employed for future clus- 
ter dynamics stud ies. Instead, the sinks in the simulations by 
IDale et al.l i|2012r ) have a mass range that goes from individ- 
ual stars to small clusters. Moreover, they only considered 
the effect of stars more massive than 20 Mq, so they did 
not investigate the effect of stars of different masses. On the 
other hand, their radiative transfer algorithm is more real- 
istic than ours, and they sample a larger parameter space. 
Thus, in general, it can be said that the two studies are 
highly complementary in nature, each one providing a dif- 
ferent perspective of the problem: specifically, they focused 
on the ability of photoinizing radiation to destroy clouds of 
different masses, while we have focused on the control of the 
SFE and the role of stars of different masses and different 
feedback prescriptions. 

Most importantly, our simulations are relevant in the 
context of studying the entire evolutionary cycle of molecu- 
lar clouds, and showing that star-forming GMCs as a whole 
can be in a global state of gravitational contraction, which is 
initiated during their pre-molecular stages, and yet comply 
with the low observed SFR and SFE in The Galaxy, as a 
consequence of the stellar feedback, but no t by maintaining 
them hovering around an equ ilibrium state (|Krumholz et al I 
120061 ; iGoldbaum et al.ll2011l ). but rather by photoevaporat- 
ing them while the collapse motions continue. 
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5.2 Limitations 

Our simulations, although one step ahead of our previous 
effort from Paper I, are still far from being all-inclusive. 
Most notably, we have neglected supernova explosions and 
magnetic fields, and moreover, our radiative transfer scheme 
is very rudimentary. We plan to improve on these issues in 
future work. Here we can speculate what should be the effect 
of these processes on our results. 

From Sec. 14. fl we see that the SFE (measured with 
respect to the maximum dense gas mass) in our sim- 
ulations, shown in Figs. 2HU with only one exception 
(Fig. [7]), is larger than ~ 5%, and reaches 30% in the 
case of the Central Cloud (Fig. [S|. This suggests that, 
in order to match the observed efficiencies of a few 
percent, a factor of a few percent further reduction in 
the SFE is still necessary. This can indeed be accom- 
plish e d by magnetic fields (e . g.. I Vazquez-Semadeni et all 
120051 ; IVazquez-Semadeni et al.ll201ll ; iNakamura fc Lilboolf ) 



and/or supernova explosions, which we plan to include in 
future studies. 

Concerning the radiative transfer (RT), our rudimen- 
tary PMRT scheme does not account for the real column 
density between an ionizing source and the test grid cell to 
be ionized, but only aims to represent it by taking the ge- 
ometric mean of the density at the source and at the test 
cell. Thus, if a dense clump lies between these two cells, 
our scheme will miss it, together with any shadowing effect 
it may have. Thus, our scheme may tend to overestimate 
the photoionized volume. We do not expect this effect to be 
dominant, since the volume covered by shadows is not large, 
but comparisons should be performed once a more thorough 
RT algorithm is implemented. 



6 SUMMARY AND CONCLUSIONS 

In this paper we have presented a numerical study of the 
entire evolutionary cycle of molecular clouds, starting from 
their formation by converging flows in the warm ISM, and 
concluding with their dispersion by the photoinization feed- 
back from the massive stars formed within them. Our study 
extended the one presented in Paper I in two main areas. 
First, we included a probabilistic scheme for star formation, 
which serendipitously allowed us to produce a stellar popu- 
lation following a realistic IMF, in turn allowing us to over- 
come a shortcoming of Paper I, namely that all stellar par- 
ticles (SPs) radiated with the same intensity, roughly corre- 
sponding to that of a ~ 10-Mq star. Second, we introduced 
our "poor man's radiative transfer scheme", PMRT, which 
allowed us to produce mass-dependent Stromgren spheres, 
and thus allowing to study the effect of stars of different 
masses in the dispersal of their parent clouds. 

We performed numerical simulations with initial con- 
ditions identical to those used in Paper I, but varying the 
feedback schemes, in order to quantify the difference be- 
tween our old and new feedback schemes (PMRT versus 
dumping all the energy in a single grid cell) and the ef- 
fect of including stelar populations of different mass ranges. 
The simulation with small-amplitude initial velocity fluctua- 
tions, SAFI (see Table[T]), due to the larger coherence of the 
converging motions in the warm gas, leads to the formation 



of a single, more massive cloud at the center of the grid, 
which we called "The Central Cloud". Instead, the simu- 
lations with large-amplitude initial fluctuations, generically 
denoted LAF, produce various less-massive clouds in the nu- 
merical box, away from the center of the simulation. From 
these we selected two, which we labeled Clouds I and 2. 

We showed that cylindrical regions of length and diam- 
eter equal to f pc were completely evacuated of dense gas 
on timescales ~ 10 Myr when a full IMF was included in 
the calculations, and in fact, the total dense gas mass in 
the numerical box is reduced by a factor ~ 5 in the SAFI 
simulation within ~ 20 Myr, and by a factor ~ 10 in run 
LAFI within ~ 15 Myr. Instead, when the most massive 
stars (M > 8 Mq) are not included in the simulations (runs 
LAFold and LAF8), the total dense gas mass in the simu- 
lations is hardly affected, although the SFE is reduced to 
levels ~ 20%. When stars up to 20 Mq are included (run 
LAF20), the total dense gas mass in the simulation is re- 
duced at a level comparable to that of run LAFI , but on 
a timescale almost twice as long. Thus, our results strongly 
suggest that the destruction of the clouds is accomplished 
by stars with masses M *S 20 Mq. 

Our simulations also show that star formation events 
can be completely terminated, and the dense gas completely 
dispersed, on scales ^ 10 pc by the photoionizing effect of 
the newly formed stars in those regions, while at larger scales 
the dense gas contents is decreased but not completely de- 
stroyed, and the SFR is analogously decreased but not termi- 
nated. This suggests that the stellar photoionizing feedback 
can locally disrupt the clouds and terminate SF, but this 
triggers new SF events at new locations in the clouds. 

Finally, our simulations show that the scenario in which 
large, dense and cold clouds begin to collapse even before 
they are mostly molecular, and continue doing so trhough 
their star-forming stages is perfectly consistent with the ob- 
served values of the SFE. 
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